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Abstract 

The different ground state energies of A^-pion and M-kaon systems for N +M < 12 are studied in 
lattice QCD. Tliese energies are then used to extract the various two- and three- body interactions 
that occur in these systems. Particular attention is paid to additional thermal states present 
in the spectrum because of the finite temporal extent. These calculations are performed using 
one ensemble of 2-1-1 flavor anisotropic lattices with a spatial lattice spacing as ~ 0.125 fm, an 
anisotropy factor £^ = ag/at = 3.5, and a spatial volume ~ (2.5 fm)'^. The quark masses 
used correspond to pion and kaon masses of ~ 383 MeV and ~ 537 MeV, respectively. 
The isospin and strangeness chemical potentials of these systems are found to be in the region 
where chiral perturbation theory and hadronic models predict a phase transition between a pion 
condensed phase and a kaon condensed phase. 
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I. INTRODUCTION 



An important goal of nuclear physics is to calculate nuclear properties and processes 
from Quantum Chromodynamics (QCD). The only known way to perform such calculations 
in an ab initio manner is using the non-perturbative formulation of QCD on a space-time 
lattice (lattice QCD). Much progress has been made in recent years in studying few hadrons 
systems, see for example the recent review of Ref. [T]. Many channels of meson- meson, 
meson-baryon and baryon-baryon interactions have been studied using the finite volume 
behavior of two-particle energy levels |2HS]. Ground state energies of three- pj and four- 
^ baryon systems have been computed for the first time in QCD and quenched QCD, re- 
spectively. To become a central part of nuclear physics, the successes of lattice QCD in the 
realm of few-hadron systems must be translated into the realm of many hadrons and the 
complexity frontier must be addressed. To this end, recent work has focused on numerical 
investigations in systems of up to twelve pions or kaons [HHH], and algorithms [12] have 
been constructed to greatly extend these studies. An important outcome of these studies 
was the extraction of the two- and three- body interactions in pion and, separately, kaon 
systems. The shifts in the A^-meson energy levels at finite volume were calculated numer- 
ically and then compared with the analytical expectations derived in Refs. [T3HI5]- These 
calculations also probed the relation between isospin(hypercharge) chemical potential and 
isospin(hypercharge) density, finding agreement with the predictions of chiral perturbation 
theory 

In this work, we study more complicated multi-meson systems involving both pions and 
kaons. As in the single species case, the expected volume dependence of the energy of an 
N-TT^, M-K^ system has been computed [18] and these calculations will enable us to extract 
the scattering lengths for the three two-body interactions, 0.^,+^,+ , a-„+K+, a,K+K+, and also to 
study the various zero momentum three-body interactions, parameterized as fjs^mnT, V3,TTnK, 
V3,ttKK, V3,kkk (restricted to the maximal isospin in each case). Multi-meson systems such 
as these exhibit a rich phase structure [Ml [19] and are phenomenologically relevant in a 
number of settings ranging from RHIC [20] to the interiors of neutron stars [21]. At low 
temperatures, the ground-state of a finite density system of, for example, pions is expected 
to undergo a transition from a Bose-Einstein condensed phase [IHl [22VI26] to a BCS-type 
superconducting phase as the isospin chemical potential is increased, and other interesting 
phenomena such as p-condensation may also occur [27H29] . Similar effects are expected in 
kaon systems. The phase structure of a mixed system of pions and kaons is a less well 
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investigated question. Studies in the NJL (Nambu-Jona-Lasinio) model, [301432] . random 
matrix theory [33], and SU(3) chiral perturbation theory (xPT) with non-zero isospin and 
hypercharge chemical potentials [17], suggest that pion and kaon condensed phases compete 
in interesting ways for different ranges of isospin and hypercharge chemical potentials. 

In the following, we focus on numerical analysis of lattice QCD realizations of these com- 
plex multi-meson systems and on the methodology needed to extract the two- and three-body 
interactions. We calculate the 90 possible two-point correlation functions involving A^-vr+s 
and M-K^s for all A^-|-M < 12. These results are then analyzed to extract the ground state 
energies of the appropriate quantum numbers, taking care to account for thermal excitations 
in the lattice volume. Finally, the resulting large set of energies is used to extract the seven 
underlying interactions discussed above. The extracted energies also allow us to study the 
chemical potential of these systems. This study is performed using a single ensemble of 2-1-1 
flavor anisotropic gauge field configurations generated by the Hadron Spectrum Collabora- 
tion [34J and so should be viewed as an exploratory study, testing methods necessary for 
analysis of complicated many body systems. Future work with different lattice spacings, vol- 
umes and quark masses will allow contact with experiment [55H5U] . in the case of two-body 
interactions, and with other lattice studies. 

II. MULTI-HADRON CORRELATION FUNCTIONS 

In order to extract the energies of the multi meson systems, we will study the behavior 
of two point correlation functions 

C^iv,M(t) = (o|(5^7r-(a^,t))'^(5^ir-(a^,t))''(7r+(0,t))'^(ir+(0,t))''|o), (2.1) 

X X 

where the pion and kaon interpolating operators are defined in terms of quark fields as 
TT^t = 7r+(x, t) = m(x, t)75(i(x, t) and = K~^{-x,t) = ?l(x, t)75(i(x, t) , respectively. This 
system corresponds to a system of A^vr+'s and Mi^+'s with total momentum, Ptot = 0. 

It is instructive to first consider the form of such correlation functions which can be 
generically expressed as 

c{t) = {o\d'\t)dmo) , (2.2) 

for some interpolating operators, 0^'\ with commensurate quantum numbers.^ In a system 



^ The creation and annihilation interpolating operators are not in general Hermitian conjugates of one 
another. 
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with infinite (Euclidean) temporal extent, application of the transfer matrix formalism shows 
that 

oo 

C{t) = J2 Zke-'^"' > (2-3) 

fc=0 

where the sum runs over all eigenstates with the quantum numbers of the operators under 
study. The Ek are the energies of eigenstates \E,,) and Zk = {0\d'\0)\Ek){Ek\d{0)\0) . At 
large Euclidean time, the correlator is dominated by the ground state energy of the system, 
C{t) —7- Zq exp{—Eot). However, in a system with a finite temporal extent, T, a natural 
choice of temporal boundary conditions for quark and gluon fields, anti-periodic and periodic 
respectively, will force C{t) to be periodic (assuming the operators under study are bosonic. 



as in Eq. (2.1)). The expected form of the correlation function can again be inferred: 

1 ^ oo 

C(t) = -tr (e-^^a't(t)a(O)) = jYl e'^^-'^^"e-'^-{Er,\d'\0)\E.^){Em\dmE^) , 

n,m=0 

(2.4) 



with Z = tT{e~^^). The contributions on the right-hand side of Eq. (2.4) can be separated 
into two distinct classes. The first class is defined by choosing \En) to be the vacuum state 
{En = 0) and correspond to the ground state and excited states of the system. These 
contributions are special as they persist in the T — t- oo (or zero temperature) limit. The 
remaining contributions exist solely due to the finite time extent and are known as thermal 
states [SI [ini HOHll]- These states vanish as T — )• oo because of the first exponential under 



the sum in Eq. (2.4). 

Thermal contributions can be illustrated with a concrete example. Consider 6*4,0 (i) which 
corresponds to a system with the quantum numbers of four pions. The correlation function 
is dominated by a sum of three contributions in the large t limit: 

C4.,oi^(t) = cosh {E.^tr) + Zse-^^^^^+^^^^^Z^ cosh ((^3^ - M^)tr) 

+ Z3 cosh {E2nT) + ■■■ 

= Zi cosh (E^Jt) + Z2 cosh ((^3^ - M^)tT) + Z3 + ■ ■ ■ , (2.5) 

where tx = t — T/2, the ellipsis denotes contributions involving excited pion systems that we 
will ignore in the current discussion, and the Z's are constant with respect to time. The first 



term in Eq. (2.5) corresponds to states where all four 7r~'s propagate forward in time {\En) = 
|0), \Em) = |47r)) and four tt+'s propagate backward in time {\En) = |47r), \Em) = |0)). The 
second term represents three vr^'s propagating forward in time and one vr^ propagating 
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backward in time {\En) = l^r), \Em) = |3vr)) as well as three vr+'s propagating backward in 
time and one tt"*" propagating forward in time {\En) = |37r), \Em) = |vr)). Finally, the last 
term arises from two tt^'s propagating forward in time and two tt^'s propagating backward 
in time = |27r), \Em) = \2n)). 

The general form for the correlation function of A^-pions and M-kaons can be straight- 
forwardly worked out and is given by the following 

M N ^ ^ 

CnM^^) ~ ^nrn^'^^ ^ COsh i^^E^^'^ "^^t) + Z ^' ^^N mod mod 2,0 + " " " ! 

m=0 n=0 

(2.6) 

where AE"^""'^"™ = {Eis[-n,M-m. — En,m), the ellipsis denotes excited state contributions, 



and the last term in Eq. (2.6) is only present when N,M are even. It is apparent that 
the number of possible terms contributing to CN,M{t) grows with and M. Locality of 
the transfer matrix guarantees that the eigen-energies En^m appearing in many places in 
different CN,M(t) are identical. Consequently, multiple correlation functions (choices of 
and M) can be used to extract the common set of eigen-energies. 

Thermal states contribute to most lattice calculations and need to be considered in pre- 
cise analyses. They are particularly prevalent in the multi-hadron context as these systems 
easily factorize into multiple color singlet states that can propagate over long distances. In 
Refs. pHTO]. the issue of thermal states was avoided by using Dirichlet boundary conditions. 
In principle, such boundary conditions introduce unknown contamination into correlation 
functions, however these works analyzed correlations significantly separated from the bound- 
ary to reduce such effects. 



III. DETAILS OF LATTICE CALCULATION 



The lattice calculations presented in this work are based on a single ensemble of 
anisotropic gauge field configurations generated by the Hadron Spectrum Collaboration |3l] . 
These configurations were generated using a stout smeared j45j gauge action and a 2+1 flavor 
clover fermion action; further details can be found in Refs. [M]. The spatial and temporal 
lattice spacings are = 0.123 fm and at = ds/C, respectively, where ^ = 3.5 is the renor- 
malized anisotropy. The lattices have a volume of 20^ x 128, corresponding in physical units 
to (2.5 fm)^ X 4.6 fm, and the light and strange quark masses are such that the pion and 
kaon have masses of ~ 383 MeV and rriK ~ 537 MeV. These configurations have been 
extensively studied by the HSC and NPLQCD collaborations. For this study, an ensemble 
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of 400 configurations were cliosen from a long stream of ~10,000 trajectories and are well 
separated in Monte-Carlo time such that autocorrelations are reduced. 

Our analysis makes use of light and strange quark propagators generated by the NPLQCD 
collaboration. For each configuration, randomly positioned APE |l6l HT] smeared sources 
were used to generate approximately 75 propagators. The EigCG inverter ^8] was used to 
perform the multiple inversions efficiently. These propagators were then APE smeared at 



the sink and combined into the correlation functions of Eq. (2.1). 



Naively, the number of contractions involved in the correlation functions of Eq. (2.1) is 



enormous, (A^ + M)\N\M\. To compute the requisite contractions, we extend the methods 
developed in Ref. [S] to the mixed species case. Following the manipulations of Ref. [U] it is 
straightforward to show that for arbitrary 12 x 12 matrices, 11, K: 



det 1 + \pU + XkK 



exp ( tr I ^ (ApH + XkK)' 



3 



oo 



-1 1 



k=l 



A;!(j-A:)! 



exp \ ti y ^ 

A^ 

1 + Ap tr n + Aa' tr ir + ^ 



tr n 



tr n- 



+ ^(^trfs: - tr fsT^ j + ^^^^ j^tr n tr - tr n/sT 



+ 



(3.1) 



with 



n(i) = J]5,(^K,t;0,0) S\{x,t-^,^), 

X 

K{t) = Y,Su{x,t;O,0) Sl{x,t;O,0), 



(3.2) 

(3.3) 



where the Sj correspond to the source and sink smeared quark propagators of flavor /. 
Additionally, one can generalize the results of the single species case to show that in the 
mixed species case, 

12 j 

■ '"^ iA^A^^-'^C,,,_fc(t), (3.4) 



1 



det (i + XpU + XkK)=^ J2Y1 

■ j=l k=0 



12 \ (3 \ y^k xi-k, 



where 



Ck,j-k(t) 



ai...Ofc/ii.../ij_S;^i...5l2-j 



e/3l...^fc!/l...!/j_fe5l...5 



nil ...UPJ^K^.I ...K, 



(3.5) 
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By expanding the exponential to a particular order in Ap and Xk-, and equating Eq. (3.1) 



and Eq. (3.4), one can identify the function Civ,j\/(^)- 

The computational complexity of these two-species contractions is significantly more than 
that of the single species cases studied previously. The reader can find an example, the 



CiT,^2,K{t) correlator, in Eq. (Bl). As in the single species cases [9l [TO], high precision 
arithmetic is required to perform these contractions correctly and this is implemented using 
the arprec and qd libraries [l9]. Explicit calculations for all + M = 13 correlation 
functions show these correlators vanish to the requisite precision and are an effective check of 
the correctness of our code. The computational cost of computing all contractions CM^uif) 
for N + M < 12 and all < t < 128 is approximately twenty minutes on a single core, 
compared to a few seconds for the single species case. While not available at the time that 
the current calculations began, the recursive constructions developed in Ref. [12] provide 
an alternate, and computationally more efficient, way of generating the contractions. The 
results of both methods agree to arbitrary precision.^ 



IV. VOLUME DEPENDENCE OF MULTI-HADRON ENERGIES 



As discussed above, the calculations of the correlators in Eq. (2.6) determine the energies 
of the mesonic systems. These can in turn be used to determine the interactions through the 
well-known results of Liischer [3l in the two body case and the results of Refs. [T3HT5| [18] 



in the many meson case where a perturbative expansion in the inverse volume is performed. 



A. Volume dependence of two particle energies 

An extraction of the scattering lengths from the single species and mixed species two 
particle systems provides a baseline reference with which to compare the results of the 
multi-particle analysis. In the center-of-mass frame, we define the interaction momentum, 
p = |p|, from the energy shift AEab = Eab — Ea — Eb = a/p^ -|- m\+^yp'^ + m\ — mA — rnB- 
AEab is determined from lattice calculations and one first solves for the momentum p [1] . 
Liischer's formula relates this momentum to the phase-shift, S{p), of two particle scattering 

^ In future calculations that extend the number of pions and kaons beyond twelve, these recursive techniques 
will be critical as the explicit code generated for all the requisite contractions would be unmanageably 
large. 
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as 



where j G and 



pcot6{p) = s[^), (4.1) 



'^(^)= E ^ttA,, (4.2) 

^ — ^ 7 — X 



is a regulated three dimensional zeta function. By expanding the left-hand side of Eq. (4.1) 
using the effective-range expansion, pcot 6{p) = — ^ + ^p"^, the scattering length, a, can be 
determined. In the absence of interactions, S possesses poles at j = 27rn/L for n G Z^. 
With interactions these poles are shifted. One then computes a as: 



(4.3) 



Using this formula with the associated one-body and two-body energies will yield precise 
determinations of the scattering lengths as the one- and two- body energies are the most 
cleanly determined and the above results only omit exponentially suppressed finite volume 
effects. 

B. Volume dependence of multi- meson energies 



In Ref. fH], the energy shift of a system A^-pions and M-kaons in a finite volume from 
the corresponding non-interacting system was calculated. The shift is given by: 

AE[n, m, L) = E{n, m^L) — n — m mx 

= AE^{n, L) + AEK{m, L) + AE^xin, m, L) , (4.4) 

with mT^K = TnT^mx / {ttIt^ + ttlk) and i G vr, if and 

47ra,- 



AEi{n,L) 



n 



rriiL^ \ 2 

— 2 



+ (M*|M+o(L-T), (4.5) 
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+ 



X' + J 



_l + 2^fn-l)(l + ^^ 



+ 2 ^ (m — 1)11 



O'-kK 



\tiL) \ \ a^K ) 

+ E E (/f ^''^^ + /f ''^ 



+ 



J=0 p=n,K ^ 

nm{n - l)r/3,^^i^(L) nm(m - l)f]3^^KK{L) 



2L6 



+ 



2L6 



(4.6) 



where the parameters a^j are related to the scattering lengths and effective ranges through 
PUB]: 



*'7r7r ^TTTT 



L3 



L3 



27ra3^r^x 
Ott/^ = a^K • (4.7) 



It is the dij parameters that will be determined in the current lattice calculations. The 
four volume dependent (but renormalization scale independent) quantities characterizing 
the momentum independent three-body interactions are defined by {y = m-^/mK)'- 

mML) = VsMf^) + ^TZ^(3v^ - 47r)log(/xL) - ^(2g[l, 1] + r[l, 1]), 



Aa 



m,-KKK{L,y) = r]3,nKK{fJ',y) " 



4at 



K 



i=0 p=TT,K AfeAfi 



EE E 

i=0 p=TT,K J\feJ\f2 



ar. 



(4.8) 



and 



M = |g(l, y), Q{y, 1), R{y, 1), R{l/y, 1/y) ^ 

U2 = {g(l, l/y), Q(y, y), R{y, y), i?(l, l/y) } , (4.9) 

with i = 7i,K. The functions Q, R, q, and r along with the coefficients fi are defined in 
Ref. [18j. The finite parts of Q{a,b) and R{a,b) are scheme dependent quantities where 
changes in the value will be compensated by changes in f]3{^); the numerical values for 
the Minimal Subtraction (MS) scheme are given in Ref. |I8|. However, the 773 are scheme 
independent and these are the quantities that will be determined during a lattice calculation. 
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Furthermore, the three-body interactions in the ttttK and ttKK cases depend on the mass 
ratio, niT^/mK- Finally, in the limit of — ?■ or M — > 0, this result simplifies to the 
previously determined A^-boson case while the limit M ^ N with tuk and all 

interactions set to be equal it simplifies to the 2A^-boson case. 



V. ANALYSIS STRATEGIES 



After calculating the above correlation functions on the ensemble of gauge configurations. 



we obtain lattice measurements of CN^uif) in Eq. (2.1). Provided that we can reliably deal 
with the presence of thermal and excited states, these allow us to extract the ground state 
energies of the A^-pion, M-kaon systems. By matching these measurements onto the analytic 



expectations of Eq. (4.6), the various two- and three-body interaction parameters can then 
be extracted. There are many ways in which such an analysis could proceed. One could 
attempt to directly fit the correlation functions in terms of overlap factors, and the scattering 
parameters and three-body interactions by inserting the explicit form of the energies from 



Eq. (4.6). However, we choose to perform this as a two stage analysis, first extracting the 



various energies by fits to the correlation functions using the model functions of Eq. (2.6) 



and then subsequently fitting these energies to the analytic forms of Eq. (4.6) to extract 
the interaction parameters. The primary advantage of this approach to performing the fits 
is increased stability, however, care must be taken to preserve the significant correlations 
between the different correlation functions which are computed on the same set of gauge 
configurations. 



A. Wavefunction Overlap 



As shown in Eq. (2.6), the Z-factors for each term are in general distinct. This makes 
for a difficult analysis since there are large numbers of linear and non-linear fit parameters 
that must be determined. In the non-interacting system, these Z-factors are simply related. 
In the single-species case, the relation would be Z^~^ = {'^)Zn such that is common 
to each term. The factorial simply counts the number of ways that ra-pions propagate 
forward and (A^ — n) propagate backward. Generalizing this to the multi-meson case leads 
to ^^m''^ = (^) (^)Z7v,A/. In an interacting system, however, this relationship does not hold. 
Nevertheless, in this work, we study pion and kaon systems with relatively weak interactions, 
and therefore we employ the ansatz above for the wavefunction factors. We stress this is a 
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crucial assumption, allowing us to subsequently fit all 90 correlators. As a check, analyses 
were also performed allowing the Z-factors to remain distinct for as many correlators as 
possible. Results from both methods were found to agree within uncertainty, but a full 
analysis using independent ^^""'^'^"'"'s proved unfeasible. 

An alternative method that can be utilized for these systems is the method of Variable 
Projection (VARPRO) |50l [51] which eliminates linear fit parameters by performing their 
minimization analytically. This approach is useful in the case of unrelated Z^~"'*^~'"'s from 



Eq. (2.4). Additionally, when many energies exist to be fitted, convergence to the minimum 
of parameter space is not always guaranteed, particularly if the minimum is flat. To aid 
convergence, Bayesian priors can be implemented so the physically reasonable region 

of parameter space is immediately tested. In an earlier stage of our analysis, the VARPRO 
method (with priors) was used for correlators with lower total particle number. The results 
were found to agree within uncertainty with the results using the ansatz above. 

As a final check, the fitted energies were substituted back into the form for the correlator 
and compared against the original data set. Consistency was achieved within uncertainty. 
Therefore, for the remainder of this work, we assume the simplified form for the wavefunction 
overlap. 



B. Thermal effects of Correlation Functions 



A particular multi-meson correlator depends not only on its ground state energy in the 
large-T limit, but also upon fewer body energies in the thermal states that were discussed 
earlier. To account for this, we began by fitting the one pion(kaon) correlator to determine 
Eifl{EQ^i) with a bootstrap method as will be explained below. Once these values were 
known, we turned to the two pion(kaon) correlator. This correlator depends on £^2,0 (-^0,2) 
and Eifi^Eo^i). Instead of refitting both energies, we used previous bootstrapped values of 
Ei o{Eq i) in the two pion(kaon) fit function. This is reasonable since the cleanest signal for 
the Ej\f^M energy will come from its respective correlator. To preserve all covariances in the 
data, the same bootstrap sample of gauge field configurations was used for all correlators. 
This sequential procedure was used to build up the analysis and ascertain all of the multi- 
meson energies. The fitted results for the energies and their statistical and systematic errors 
can be found in Tables |T] and Ull 
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C. Statistical Analysis 



The correlation functions span the time interval [0, 127] and a suitable subinterval, where 
the fit formula is applicable, must be chosen. As a first step, the (periodic) correlation 
function data was reflected around its midpoint in order to improve statistics. Effective mass 
plots for all correlation functions were constructed, allowing initial estimates to be made of 
where (in time) contamination from excited states ceases. These plots revealed that excited 
state contributions are suppressed by t ~ 20. Therefore, in an effort to be as general as 
possible, fits were run on all time intervals, (tmin,^min + '^t), such that tmin & [25,33] and 
At G [10,20]. Using intervals beyond these values showed no improvement in fit results. 

In this work, we utilized the bootstrap method of statistical resampling to estimate 
uncertainties. To construct the bootstrap ensembles, we first averaged over correlators from 
all sources on each gauge configurations, for each N,M. Let Ck{ti) be a given correlation 
function for some fixed number of pions and kaons, N and M, calculated on a given gauge 
field configuration k at time tj. Then we denote this set as {Ci(tj), C2(ti), ■ ■ ■ , C_A4(ti)}, with 
Ai = 400. The data covariance matrix is defined from the full data set as: 



1 



M 



^(^-^^■) = A.M. 1^ E (^'^(^^) - ^(^^)) - Cit,)) , (5.1) 



M{M-l) 



k=l 



where C{U) = (l/M) ZZi Ck{U). 

Next, we randomly sampled configurations from the full data set and constructed a 
bootstrap sample, -Bi(ti) of Ai elements, for each time ti. This sample was averaged 
to yield -Bi(tj) and this entire process repeated V times where V = 450, resulting in: 
{Bi{ti), B2{ti), ■ • • , B-p(ti)}. For each of the Bj, a x] function was determined as: 

^? = E {B,{t)-y{t;EN,M)) C-\t,t') (^B,{t') - y{t'; E^,,,)) , (5.2) 
t,t' 



where y{t; E^^m) represents the fit model in Eq. (2.6). Each of the Xj functions were mini- 
mized and results recorded. 

In order to assess the reliability of the fits, we define the goodness-of-fit, Qj = 
2d/2Y(d/2) I(x^ ) '^^i (Xj)'^^^~^6~^?^^, where d is the number of degrees-of-freedom, for each 
Xj- Once a best fit energy is selected on each bootstrap sample, we determine a weighted 
average for the ensemble, 

1 ^ 

{E^,m) = y^Y.Q^^N,M, (5.3) 
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with weights Qj, and a weighted standard deviation, 



a 



EQ.(^S!a/-(^^,m))'. (5.4) 



Above, Wi = E;=i and W2 = E,=i Q'- 

To assign a systematic uncertainty for the extracted energies, we return to the fits per- 
formed to the correlators. This is a necessary component of the analysis since there exists 
arbitrariness in the choice of which tmin and which At were selected for the correlation 
function fits as explained above. The systematic error on the best fit energy is defined as 

1 



E^Mi^U) - ^^,Af(At_) , (5.5) 



where -E'Ar,Af(At+,_) are the energies extracted from CN,M(t) using the time intervals [tmm + 

1; tjYiax ~l~ 1] and [tmin 1; tjnax 1]. 

D. Two and Three Body Parameters 

The second stage of our analysis now focuses on the extraction of the scattering lengths 



and three-body coefficients from the measured energies using Eqs. (4.5) and (4.6). These 
determinations rely on a second bootstrap analysis involving a resampling of the extracted 
energies. The bootstrapping procedure for a specific correlation function yielded V ener- 
gies, and these formed the bootstrap samples for the extraction of the two and three-body 
parameters. 

Once the best fit multi-meson energies were known, a very similar procedure used for 
the analysis of the correlators was used to find the a's and f/'s. Since a bootstrap ensemble 
exists for every best fit energy value, we created an energy sample, £a such that a G [Ij'P]. 
This sample carries an additional vector index that labels the energies within the vector. In 
the case of single species pion energies (the kaon case is identical), an energy vector initially 
composed of Sa = {-^2*0* 5 ^30 ^ ^ao} *as used to fit to Otttt and f]^,Tvmv We included another 
energy and refitted the interaction parameters and repeated this until all the energies were 
exhausted. In the multi-species case, a base set of {-^2^0^ ' ' ' ' -^12 05 ^0^21 " 1 -^o"i2} along 
with ten randomly selected multi-species energies was created and fits performed for all 
seven hadronic parameters. This first set thus made use of 34 different energies. This set 
was enlarged by one, the parameters were refitted, and the process repeated until all ninety 
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energies were used. The energy covariance matrix used in these fits, is defined according to: 

1 ^ 

= E (^-.^ - (^')) (^"-^- - (^^■)) ' (5.6) 



a=l 



such that {£i) = (l/P) J2a=i^a,i, and the energy on each bootstrap is defined as: 



x: 



where f{ci,ff) is shorthand notation for the fit functions in Eqs. (4.4)-(4.6). 

The systematic errors assigned to the a's and f^'s are more comphcated than those of the 
energies. Given a particular energy set of M energies that are used to make a determination 
of a's and f^'s there are 3-^ different combinations of the intervals that must be fit in order to 
completely propagate the systematic uncertainties of the energies to those of the interaction 
parameters (it is 3-^ because there is a [tmin^tmax] for each best fit energy as well as its sys- 
tematic counterparts corresponding to the shifted time interval in the forward and backward 
direction). Even in the single species case, when M = 10, there are already ~ 6 x 10^ com- 
binations. For the multi-species case, it is too costly to fit all these permutations. Rather, 
we only fit O{10^) randomly chosen permutations and take the difference of the mean of 
this set from the best fit a and f] as the systematic error. From fitting all permutations in 
the single species case, up to A/" = 9, it was seen the systematic error stabilized well before 
the total number of combinations was computed and we assume this is also the case for the 
two-species case. 

VI. RESULTS 
A. Energies 

Using the methods discussed above, we extracted the energies of the mixed and pure 
species system, from all ninety correlators. The final extracted values are shown in Tables 
|I] and |Tl] below, along with their associated fit ranges. These energies are shown in a three- 
dimensional plot along with their respective uncertainties in Fig. [TJ 

The fits become progressively more difficult as the number of mesons grows because of 
the increasing thermal contamination. This is directly refiected in the quality of the fits 
decreasing for large meson number in both the pure species and mixed species case. Fits to 



example correlators are shown in Fig. VIA 
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FIG. 1: Energy of multi meson states. Uncertainties shown are result from combining statistical 
and systematic uncertainties in quadrature. 

B. Interactions 



The extractions of interaction parameters from mixed meson energies were performed to 
yield the three scattering lengths and four three-body coefficients. This work builds upon 
the studies of [HHIS| and presents the first measurements of fjs^mrK, and f/Trxx since these 
parameters can only be measured within the framework of the mixed-meson system. 

The most straightforward determination of the scattering lengths is given by using the 



eigenvalue relation from Eq. (4.1). Using this, we find 



m^a^^ = 0.225 ± 0.001 ± 0.023, 
TJiKaKK = 0.4465 ± 0.0006 ± 0.0266, 
m^KanK = 0.1560 ± 0.0004 ± 0.0095. 



(6.1) 
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FIG. 2: Plots of the log of the fitted correlation function (red) and those based on the full data set 
of gauge configurations (blue) for their respective fit intervals for representative N and M. The 
red envelope denotes the uncertainty in the fitted correlator propagated from the uncertainty of 
the energies. 



The three-body coefficients can only be determined within the framework of Eqs. (4.4), (4.5) 



and (4.6). We also use this same analysis to provide a check on the above results. Given 
that our analysis provides multiple determinations of the interactions parameters for varying 
numbers of combinations of energies used in the fits, these must be combined in some way 
to obtain the final values. Since each separate extraction can be viewed as a somewhat 
independent measurement, the final value given is taken to be the mean from the set of 
all extractions. The final uncertainties on the extractions are combinations of statistical 
uncertainties, systematic uncertainties obtained from variation of the fitting windows as 



discussed in Sec. |VD| and a second systematic uncertainty determined from the standard 
deviation of the full set of extractions, combined in quadrature. The systematics are the 
largest source of uncertainties in the results. The individual extractions of the various 
parameters and the final extractions are shown in Figs. IMTI The error bars shown combine 
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the statistical, and systematic uncertainties as discussed in Sec. |V] in quadrature. The 
shaded regions with thin borders denote the final results and their uncertainties. For the 
mixed species extractions, the second shaded band with thick borders, denotes the range of 
uncertainty in the quoted values from the single species analysis. These are shown together 
so the reader can see the overlap region between both sets of results. The poorest behavior 
originates from fj^^Trnn where the mixed species results drift away from the pure species one. 
The final values of the interaction parameters for the single-species case are: 

mxaKK = 0.444 ± 0.011, m^a^^ = 0.224 ± 0.031, 
mKr]3,KKKfK = 0-11 ± 0.28, m^r/a,.../,' = 1.81 ± 0.52, (6.2) 



whereas for the multi-species case we find: 



rriKaKK = 0.461 ± 0.010 



m^a^^ = 0.271 ± 0.021, 



fTT-nKClnK = 0.166 ± 0.016, 

mKm,KKKfK = -0.08 ±0.12, 



"^7r?73,7r7r7r/^ = 0.68 ± 0.33, 

r/3,.ifi^/'xi^ = 0.22 ±0.17, 



m^r ± 2mK 

^;^,3,../l..^ 0.45 ±0.26. (6.3) 

For the two-body parameters, the perturbative analysis (for both pure species and mixed 
species) is in agreement with the nonperturbative values. In the above results, only the 
final uncertainty, including all statistical and systematic contributions, is given as discussed 
above. A full list of our extractions can be found in Table IIIIl 

In previous studies, NPLQCD [54j, CP-PACS |55], and ETMC |56] have measured the 
pion scattering lengths in the isospin, 1 = 2 channel. The NPLQCD determination of the 
7r+7r+ scattering length at m^r = 350 MeV, yielded |m7ra7r7r| = 0.2061(49)(17)(20) whereas an 
analysis from multi-pion correlators [H |9] yielded an extraction of \mTrO.Tnr\ = 0.2058(45)^82 
and mT,fi2,^T,T,T, = 1.02(08)^22 (i^i these results, the first uncertainty is statistical and the 
other uncertainties are systematics as discussed in the original references. In the ETMC 
analysis of 7r+7r+ scattering, |m,rfl7r7r| = 0.252(22) (13) is measured at m,r ~ 391 MeV. 
Fig. Isl shows the dimensionless combination |m7ra^7r| from the current work in compari- 



son to the determinations by other groups at a similar pion mass. However, these results 
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FIG. 3: Calculation of scattering lengths and three-body coefficients for kaons from Eq. (4.4). Un- 
certainties per extraction result from combining statistical and systematic uncertainties in quadra- 
ture. The shaded band defines the standard deviation of the mean of all extractions. The red point 
is the nonperturbative Liischer result. 



are at non-zero lattice spacing and correspond to different discretizations, so agreement is 
not necessary. The / = 1, K^K^ scattering length was also determined by the NPLQCD 
collaboration in Ref. [B7J. An analysis of the two-point kaon correlator yielded a value of 
\tt^kO'Kk\ = 0.497(10) (22), again at ~ 350 MeV. Analysis of multi-kaon correlators [TU] . 
led to \mKaKK\ = 0.503(11)(19) and mKfjz,KKKfK — ~0-l(2)(5) where the uncertainties 
are statistical and systematic respectively. The tcK scattering length has been investigated 
in quenched |SB] and full QCD [SS] and the unquenched determination at ~ 350 MeV 
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FIG. 4: The meaning of the points is identical to that of Fig. |3] only we show pions rather than 
kaons. 

is rriT^KO'TTK = 0.155(40). Hence, it is clear the current results are generally consistent with 
other groups' extractions. The mixed species three-body parameters are novel results and 
are found to be of natural size and positive. 



C. Isospin and Hypercharge Chemical Potentials 

As we have determined the dependence of the energy of the mixed meson systems on the 
number of pions and kaons, we can construct the isospin and hypercharge (or strangeness) 
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FIG. 5: Calculation of hadronic parameters for mixed- mesons. Extractions are shown as a function 
of the number of energies that were used in the extraction. The number of energies used ranged 



from 34 to 90 as discussed in Sec. VD Uncertainties for a given extraction result from combining 
statistical and systematic uncertainties in quadrature. The shaded band with thick borders denotes 
the standard deviation of the mean of all extractions in the multi-species case while the shaded band 
with thin borders denotes the standard deviation of the mean of all extractions in the single-species 
case. The red point at the left-most end is the nonperturbative Liischer result. 



chemical potentials using finite differences following Refs. [HI [ID] where systems of pions 
and kaons were investigated separately. In Refs. [HI [ID], remarkable agreement was found 
between the numerical results and the leading order (LO) xPT prediction [60j of the relation 
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FIG. 6: Scattering parameters and three-body interactions are shown. The meaning of the points 
and regions is the same as in Fig. [5} 



between the isospin(hypercharge) density and chemical potential 



m 



fx 



4 /' 



(6.4) 



with j G {7c,K). The situation here is more complicated since there are finite differences 
acting in various non-orthogonal directions; the differences between E^^^^^^ and E^^_i^^j^ 
determine /i/ while linear combinations of E^^^Nk^ ^Nt,-i,Nki -EAr^,Ar^-i determine ns^- 



fiK+ — fJ'S + with — En„,Nk ~ Ej\[^,Nk-i and /i/ — En„,Nk " ^JV^-i,Afj< 
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FIG. 7: Scattering parameters and three-body interactions are shown. The meaning of the points 
and regions is the same as in Fig. [5| 

One goal of this analysis is to see where on the fis vs. phase diagram [17j, the states 
created in the lattice calculation lie. 

In Ref. [T7], leading order SU(3) xPT is used to predict three distinct phases for nonzero 
isospin and hypercharge chemical potential. The first is the normal phase where the ground 
state has a net particle number of zero. The other two phases are the pion-condensed and 
kaon-condensed phases. The transition between the kaon-condensed phase and the pion- 
condensed phase is predicted to be a first order phase transition, separated by the line 
fxs = ("""^^ + a/ {rnl — + Amj^^'j ) /2/x/, while the transition from the normal phase 
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0.18 0.20 0.22 0.24 0.26 0.28 0.30 

FIG. 8: The values of |m7ra,r7r| obtained by different groups with pion masses, ~ 350 MeV are 
shown. From bottom to top, the data are from NPLQCD NPLQCD [SlIS], ETMC [56], the 
present work's single species value, and the present work's multi-species value, respectively. Note 
that these calculations are at non-zero lattice spacing and use different discretizations so complete 
agreement is not expected. 

to either condensed phase is expected to be of second order and are defined by the lines 
fis = ~ yU//2 and fij = m^. These predictions assume zero temperature and are likely 
softened by the non-zero temperature at which the lattice calculation is performed [62j . 

In Fig.|9| both the lattice calculations of (/x/, /i^) and the xPT phase boundaries are shown 
(dashed lines). Data points corresponding to higher numbers of particle states are shown in a 
orange/reddish color, while lower numbers are given in a blueish/greenish color. Points with 
large uncertainties are excluded from this figure for clarity (the omitted data correspond to 
the highest particle numbers). It is striking that the calculated chemical potentials mostly 
lie near the first-order phase transition line predicted by xPT. Further calculations with 
larger numbers of pions and kaons will be enlightening, but more complex probes of these 
systems may be needed to fully understand the states that have been produced. 

^ An AdS/QCD based model [ST] finds these transitions to be of first order. 
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FIG. 9: us ys. iJ,[. The data points corresponds to the lattice calculations, colored by the energy 
of the system from low (blue) to high (orange). The dashed lines are predictions of xPT. The 
lower-left region is the normal phase, the right-hand region is the pion- condensed phase, will the 
upper portion is the kaon- condensed phase. 



VII. CONCLUSIONS 



In this work, we have numerically studied complex systems of mesons of two distinct 
flavors, like-charged pions and kaons, and used them to extract information about the two- 
and three- body interactions amongst pions and kaons. Where known, the interactions were 
found to be consistent with previous calculations, however, two mixed-species three-body 
interactions were determined for the flrst time. Additionally, the isospin and strangeness 
chemical potentials and phase structure of the system have been investigated, with the 
systems preferring to probe a region in the (/i/,/Xs) plane where xPT predicts a flrst order 
phase transition. 

A major aim of this work was to investigate technical issues that arise in the analysis of 
complex multi-hadron systems. Accounting for the thermal states that proliferate in such 
systems, which easily factorize into distinct color singlet states, proved challenging and future 
calculations should avoid this by using larger temporal extents. Additionally, a number of 
techniques to perform coupled flts to the C(IOO) correlators studied were investigated and 
found to be beneflcial in the analysis. 
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In the future, calculations probing larger meson numbers will allow further investigations 
of the phase structure of these interesting QCD systems. To understand the structure of 
the condensed systems created in the current and future calculations, more complicated 
observables that access transport properties may be needed; investigations in this direction 
are under consideration. 
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Appendix A: Tables 



TABLE I: Extraction of single-species meson energies (in lattice units). The first set of uncertainties 
are statistical and the second set are systematic uncertainties associated with shifting the correlator 
fitting window by At = ±1. The remaining column shows the fit range chosen. 



AT Ayf 

Jv, M 


/ C^X _l_ _ _l_ _ 

{-tj) ± CTstat ± CTsys 


(i- f \ 
\'^min 1 '•max J 


1,0 


0.06936 ± 0.00021 ± 0.00002 


(27,37) 


2,0 


0.14067 ± 0.00049 ± 0.00007 


(27,37) 


3,0 


0.2143 ± 0.0010 ± 0.0004 


(27,37) 


4,0 


0.291 ± 0.002 ± 0.001 


(27,37) 


5,0 


0.378 ± 0.003 ± 0.003 


(31,41) 


6,0 


0.470 ±0.005 ±0.004 


(31,41) 


7,0 


0.547 ±0.005 ±0.003 


(26,36) 


8,0 


0.654 ±0.007 ±0.003 


(31,41) 

\" 7 / 


9,0 


0.7560 ± 0.0090 ± 0.0007 


(31,41) 

\" 7 / 


10,0 


0.868 ±0.014 ±0.002 


(31,41) 


11,0 


0.983 ±0.015 ±0.003 


(32,42) 


12,0 


1.101 ±0.046 ±0.009 


(25,35) 


0,1 


0.09727 ± 0.00015 ± 0.00001 


(27,37) 


0, 2 


u.iyoDo ± u.uuuoi ± u.uuuu( 


(27,37) 


0,3 


0.2982 ± 0.0005 ± 0.0004 


(27,37) 


0,4 


0.4021 ±0.0008 ±0.0011 


(27,37) 


0,5 


0.5083 ±0.0011 ±0.0008 


(27,37) 


0,6 


0.617 ±0.002 ±0.003 


(27,37) 


0,7 


0.728 ±0.002 ±0.007 


(28,39) 


0,8 


0.842 ± 0.003 ± 0.004 


(28,39) 


0,9 


0.959 ±0.005 ±0.011 


(28,39) 


0,10 


1.101 ±0.036 ±0.009 


(32,49) 


0,11 


1.23 ±0.03 ±0.01 


(27,41) 


0,12 


1.384 ±0.076 ±0.001 


(33,52) 



TABLE IL Extraction of multi-species meson energies (in 
lattice units). The columns are the same as in Table |lj 



1,1 
1,2 
1,3 
1,4 
1,5 
1,6 
1,7 
1,8 



{E) ± (Tstat ± cr. 



0.1687 ±0.0004 
0.2703 ±0.0006 
0.3743 ± 0.0009 
0.481 ± 0.001 ± 
0.590 ± 0.002 ± 
0.7024 ±0.0032 
0.8171 ±0.0043 
0.938 ± 0.005 ± 



± 0.0002 
± 0.0006 
±0.0013 
0.002 
0.004 
± 0.0005 
± 0.0003 
0.005 



It 



tr 



mm, "-max 



J 



27,37 
27,37 
27,37 
27,37 
27,37 
28,39 
28,39 
27,37 



Continued on next page 
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TAB 


LiE II continued from previous page 


N,M 


{E) zb CJstat ± CTsys 


(^min) ^max) 




1.07 ± 0.01 ± 0.12 


31,49; 




1 in 


1.21 ± 0.02 ± 0.19 


27,4i; 




1 11 

1 , i 1 


1.36 ± 0.08 ± 0.05 


'33,49; 




9 1 


0.2423 ± 0.0008 ± 0.0005 


27,37 




9 9 


0.347 lb 0.001 ± 0.001 


,27,37, 






0.453 ± 0.002 ± 0.002 


'27,37; 




9 4 


0.563 ±0.002 ±0.003 


.27,37^ 




9 ^ 
z, o 


0.6770 ± 0.0037 ± 0.0007 


'28,39; 




9 fi 

z, u 


0.793 ± 0.006 ± 0.003 


^28,39, 




9 7 


0.914 ± 0.006 ± 0.007 


'27.37; 




9 8 

Z, o 


1.049 ± 0.006 ± 0.057 


.31,49; 




9 Q 

Z, y 


1.19 ± 0.01 ± 0.04 


'33,50; 




9 1 n 

Z, iU 


1.33 ± 0.04 ± 0.01 


,33,49, 




^ 1 


0.319 ± 0.001 ± 0.001 


'27,37; 




9 

o, z 


0.426 ± 0.002 ± 0.002 


27,37 






0.537 ± 0.003 ± 0.002 


,27,37, 




O, 4: 


0.652 ± 0.005 ± 0.002 


'28,39; 






0.769 ±0.007 ±0.004 


^28,39, 






0.889 ± 0.008 ± 0.003 


'28,38; 




7 


1.02 ± 0.01 ± 0.10 


^27,38, 




o, o 


1.16 ±0.03 ±0.01 


'33,48; 






1.305 ± 0.040 ± 0.006 


,33,49, 




A 1 


0.400 ± 0.003 ± 0.001 


'27,37; 




zL 9 


0.512 ± 0.004 ± 0.001 


^27,37, 






0.626 ±0.004 ±0.002 


'28,39; 






0.745 ± 0.007 ± 0.005 


28,39 




A ^ 


0.868 ± 0.011 ± 0.008 


,28,39, 




A R 
4t, D 


0.99 ± 0.01 ± 0.05 


'27,38; 




A 7 


1.132 ± 0.017 ± 0.009 


,33,47, 




A 8 


1.280 ±0.03 ±0.001 


'33,49; 




^ 1 


0.492 ± 0.004 ± 0.004 


,31,41; 




^ 9 
0, Z 


0.610 ± 0.006 ± 0.004 


'3i,4i; 




Oj O 


0.7313 ± 0.0061 ± 0.0007 


31,42 




A 


0.856 ± 0.010 ± 0.001 


31,42 




0, 


0.982 ± 0.01 ± 0.007 


,31,43, 




0, D 


1.113 ± 0.023 ± 0.002 


'33,46; 




^ 7 


1.250 ± 0.028 ± 0.003 


33,47 




fi 1 


0.590 ± 0.007 ± 0.003 


,31,41, 




R 9 
0, z 


0.697 ± 0.006 ± 0.007 


'27,37; 




u, o 


0.837 ±0.009 ±0.002 


,31,42; 




A 


0.968 ± 0.012 ± 0.002 


'3i,4i; 




U, u 


1.103 ± 0.025 ± 0.006 


,33,43; 




Dj D 


1.2349 ± 0.0244 ± 0.0005 


'33,47; 




7 1 


0.664 ± 0.006 ± 0.007 


^26,36, 




7 9 


0.8098 ± 0.009 ± 0.0003 


'33,43; 




7 


0.949 ± 0.011 ± 0.003 


,31,41; 




7 A 


1.091 ± 0.016 ± 0.002 


'32,42; 




7 ^ 


1.19 ± 0.05 ± 0.02 


28,38 




8 1 


0.770 ± O.OOb ± 0.009 


^27,37, 




8 9 

o, Z 


0.914 ±0.010 ±0.003 


'3i,4i; 




8, 3 


1.071 ±0.015 ±0.002 


^32,42; 




4 


1.160 ±0.051 ±0.027 


'28,38; 




9^1 


0.8757 ± 0.0084 ± 0.0002 


^31,42; 




9,2 


1.034 ±0.014 ±0.002 


'32,42; 




9,3 


1.13 ±0.05 ±0.03 


^28,38; 




10,1 


0.988 ±0.011 ±0.002 


'32,42; 




10,2 


1.11 ±0.04 ±0.03 


^28,38; 






Continued on next 


page 



TABLE II continued from previous page 



N,M\ {E)±a,tat±CT: 



sys 



I* 



mm; Tnax 



111, 1 II 1.105 ±0.014 ±0.005 



(32,42) 



Letting: 

mKf}3,KfK: F5 = m^773,^/^, 

m^rriKm^^KKfnKK p _ '^■KmKfiz,-K'KKft'KK 
+ 2mK ' 2m^ + rriK 



TABLE III: Results for scattering and three-body parameters 
in lattice units. The uncertainties shown are statistical and 
systematic, respectively. 



JNum 


Fi 


Fo 






F. 




F(, 


F7 


'\A 
0^ 


0.451 








0.22(1) 









141(9)(4) 


0.046( 


21 




1.4(16)(3) 


0.43(10 






0.83( 


14 




00 


0.452 


;1S 






0.22(1) 









1422(89)(7) 


0.056( 


21 


IS 


1.39(165)(8) 


0.43(10 


il^i 




0.83( 


14 




36 


0.452 


10 


) 4; 




0.22 1 









142 


:8)(2) 


0.056(22)(0) 


1.3(16)(1) 


0.4(1)(- 




0.82(14)(2) 


37 


0.455 


11 


) 




0.23 1 









133 


11 






0.1 








1.2 


'16 9 


0.57 11 






1. 15 


)(6 




38 


0.455 


^11 






0.245(14)(2) 





135 


vie 






0.1 


1 






1-(1)(2) 




0.55 


:i2 






1.(15 


6 


39 


0.453 


'11 






0.24(12)(9) 





131 






0.1 
















0.56 


10 












40 


0.451 


^11 






0.242(12)(3) 





133 






0.1 


1! 


11! 




1:11 


11^ 


I 




0.54(10 






1:11! 


111 




41 


0.458 


'ii)(5: 




0.26(1) 









151 








0.4 
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TABLE III — continued from previous page 
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Appendix B: Example Multi-Meson Correlation Function 



C,^Mt) = :^tr[K]hr[Ur - ^trli^^jtr [K'] tr[n]^ + ^tr [K^] tr[n]^ - ^tr[K]MU]hr [U'] 

+ ^tr[X]tr [K^] tT[U]Hr [H^] - ^tr [K^] tr[Il]Hr [if] + ^tr[X]=^tr [n^]' 

-^tr[K]tr [K^] tr [n^]' + ^tr [K^] tr [H^]' + ^tr[i^]3tr[n]tr [U^] 

--tT[K]tT [K^] tr[n]tr [n^] + -tr [K^] tr[n]tr [n^] - ^tT[K]hT [U^] 
6 9 24 

+ -tr[X]tr [X^l tr [n^l - — tr [K^] tr [n^l - —tT[K]HT[U]hT[KU] 
8 12 12 

+— tr [K^] tr[n]3tr[ii'n] + ^trf/sTj^trpjtr [U^] tT[KU] - hi [K'^] tr[n]tr [U^] tT[KU] 

-^tr[X]2tr [n=^] tr[Xn] + ^tr [K^] tr [n^] tr[Xn] + ^tr[X]tr[n]2tr[Xn]2 

-^ti[K]tT [n^] tT[KU]^ - ^tT[U]tT[KUf + ^tr[/s:]2tr[n]V [KU^] 

-^tr [X^] tr[n]2tr [XH^] - ^tr[K]hr [U^] tr [XH^] + ^tr [K^] tr [n^] tr [KU^] 

-tr[/s:]tr[n]tr[/s:n]tr [KU'^] + -tr[KU]hr [KU'^] + -tr[/s:]tr [/sTH^]^ 

--tr[X]^tr[n]tr [XH^] + -tr [K'^] tr[n]tr [i^n^^] + tr[X]tr[Xn]tr [KU^] 
2 2 

+ ^tr[ii']2tr [ii^n^] - ^tr [K^] tr [i^H^] + -tr[K]tr[U]hr [K^U] 
2 2 6 

--tr[ir]tr[n]tr [n^] tr [K^U] + -tr[ir]tr [n^] tr [K^U] - -tr[n]2tr[irn]tr [K^U] 
2 3 2 

+ -tr [n^] tr[i^n]tr [K^U] + tr[n]tr [KU'^] tr [isT^n] _ tr [KU^] tr [i^^n] 
2 

-^tr[K]tr[n]2tr [K^U^] + ^tr[/s:]tr [n^] tr [K^U^] + -tr[n]tr[irn]tr [iT^^tf] 
4 4 2 

-^tr [Xn^] tr [X^n^j + 2tr[X]tr[n]tr [K'^U^] - 2tr[Xn]tr [K'^U^] 

--tr[ir]tr \K'^Vi^\ - -tr[n]3tr [/sT^n] + itr[n]tr [H^] tr [/sT^n] 

2 6 2 

--tr [n=^] tr [x^n] +tr[n]V [x^^n^] - tr [n^] tr [x^n^] - ^tr[n]tr [x=^n=^] 

3 3 

+5tr [X^n^] (Bl) 
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